# encoding:utf-8
import cv2 as cv
import copy
import numpy as np


# 打开图像
filename = r'D:\AI\image\rice.jpg'
image = cv.imread(filename)
#转化为灰度，threshold只支持灰度图像，所以先转换
gray = cv.cvtColor(image, cv.COLOR_BGR2GRAY)
# 大津算法灰度阈值化
_, bw = cv.threshold(gray, 0, 0xff, cv.THRESH_OTSU)

# 画出灰度直方图
#plt.hist(gray.ravel(), 256, [0, 256])
#plt.show()
#形态学处理，去除噪声
element = cv.getStructuringElement(cv.MORPH_CROSS, (3, 3))
bw = cv.morphologyEx(bw, cv.MORPH_OPEN, element)
#拷贝一份图像
#以下是图像分割
seg = copy.deepcopy(bw)
# 计算轮廓
bin,cnts, hier = cv.findContours(seg, cv.RETR_EXTERNAL, cv.CHAIN_APPROX_SIMPLE)
rice_area = []
rice_w = []
rice_h = []
count = 0
# 遍历所有区域，并去除面积过小的
for i in range(len(cnts), 0, -1):
    c = cnts[i-1]
    area = cv.contourArea(c)#计算面积，滤除面积小于10的分割结果：可能是噪声
    if area < 10:
        continue
    rice_area.append(area)
    count = count + 1
    print("blob", i, " : ", area)

    # 区域画框并标记
    x, y, w, h = cv.boundingRect(c)#在原始图像上画出包围矩形，并给出矩形标号
    rice_w.append(w)
    rice_h.append(h)
    cv.rectangle(image, (x, y), (x+w, y+h), (0, 0, 0xff), 1)
    cv.putText(image, str(count), (x, y), cv.FONT_HERSHEY_PLAIN, 0.5, (0, 0xff, 0))

print("米粒数量： ", count)
cv.imshow("source image", image)
cv.imshow("otsu image", bw)

rice_area_mean = np.mean(rice_area)
rice_w_mean = np.mean(rice_w)
rice_h_mean = np.mean(rice_h)

rice_area_var = np.var(rice_area)
rice_w_var = np.var(rice_w)
rice_h_var = np.var(rice_h)

rice_area_std = np.std(rice_area)

num = 0
for i in rice_area:
    if i<(3* rice_area_std):
        num = num +1
print("米粒面积均值",rice_area_mean)
print("米粒宽度均值",rice_w_mean)
print("米粒高度均值",rice_h_mean)

print("米粒面积方差",rice_area_var)
print("米粒宽度方差",rice_w_var)
print("米粒高度方差",rice_h_var)
print("米粒面积标准差",rice_area_std)

print("落在3sisgma范围内",num)
cv.waitKey()
cv.destroyAllWindows()